State and parameter estimation from exact partial state observation in stochastic reaction networks.

We consider chemical reaction networks modeled by a discrete state and continuous in time Markov process for the vector copy number of the species and provide a novel particle filter method for state and parameter estimation based on exact observation of some of the species in continuous time. The conditional probability distribution of the unobserved states is shown to satisfy a system of differential equations with jumps. We provide a method of simulating a process that is a proxy for the vector copy number of the unobserved species along with a weight. The resulting weighted Monte Carlo simulation is then used to compute the conditional probability distribution of the unobserved species. We also show how our algorithm can be adapted for a Bayesian estimation of parameters and for the estimation of a past state value based on observations up to a future time.

Stochastic Dynamics of Eukaryotic Flagellar Growth.

We study the dynamics of flagellar growth in eukaryotes where intraflagellar transporters (IFT) play a crucial role. First we investigate a stochastic version of the original balance point model where a constant number of IFT particles move up and down the flagellum. The detailed model is a discrete event vector-valued Markov process occurring in continuous time. First the detailed stochastic model is compared and contrasted with a simple scalar ordinary differential equation (ODE) model of flagellar growth. Numerical simulations reveal that the steady-state mean value of the stochastic model is well approximated by the ODE model. Then we derive a scalar stochastic differential equation (SDE) as a first approximation and obtain a "small noise" approximation showing flagellar length to be Gaussian with mean and variance governed by simple ODEs. The accuracy of the small noise model is compared favorably with the numerical simulation results of the detailed model. Secondly, we derive a revised SDE for flagellar length following the revised balance point model proposed in 2009 in which IFT particles move in trains instead of in isolation. Small noise approximation of the revised SDE yields the same approximate Gaussian distribution for the flagellar length as the SDE corresponding to the original balance point model.

Predicting Drug Release From Degradable Hydrogels Using Fluorescence Correlation Spectroscopy and Mathematical Modeling.

Predicting release from degradable hydrogels is challenging but highly valuable in a multitude of applications such as drug delivery and tissue engineering. In this study, we developed a simple mathematical and computational model that accounts for time-varying diffusivity and geometry to predict solute release profiles from degradable hydrogels. Our approach was to use time snapshots of diffusivity and hydrogel geometry data measured experimentally as inputs to a computational model which predicts release profile. We used two model proteins of varying molecular weights: bovine serum albumin (BSA; 66 kDa) and immunoglobulin G (IgG; 150 kDa). We used fluorescence correlation spectroscopy (FCS) to determine protein diffusivity as a function of hydrogel degradation. We tracked changes in gel geometry over the same time period. Curve fits to the diffusivity and geometry data were used as inputs to the computational model to predict the protein release profiles from the degradable hydrogels. We validated the model using conventional bulk release experiments. Because we approached the hydrogel as a black box, the model is particularly valuable for hydrogel systems whose degradation mechanisms are not known or cannot be accurately modeled.

Homogenization Theory for the Prediction of Obstructed Solute Diffusivity in Macromolecular Solutions.

The study of diffusion in macromolecular solutions is important in many biomedical applications such as separations, drug delivery, and cell encapsulation, and key for many biological processes such as protein assembly and interstitial transport. Not surprisingly, multiple models for the a-priori prediction of diffusion in macromolecular environments have been proposed. However, most models include parameters that are not readily measurable, are specific to the polymer-solute-solvent system, or are fitted and do not have a physical meaning. Here, for the first time, we develop a homogenization theory framework for the prediction of effective solute diffusivity in macromolecular environments based on physical parameters that are easily measurable and not specific to the macromolecule-solute-solvent system. Homogenization theory is useful for situations where knowledge of fine-scale parameters is used to predict bulk system behavior. As a first approximation, we focus on a model where the solute is subjected to obstructed diffusion via stationary spherical obstacles. We find that the homogenization theory results agree well with computationally more expensive Monte Carlo simulations. Moreover, the homogenization theory agrees with effective diffusivities of a solute in dilute and semi-dilute polymer solutions measured using fluorescence correlation spectroscopy. Lastly, we provide a mathematical formula for the effective diffusivity in terms of a non-dimensional and easily measurable geometric system parameter.

SPSens: a software package for stochastic parameter sensitivity analysis of biochemical reaction networks.

SUMMARY: SPSens is a software package for the efficient computation of stochastic parameter sensitivities of biochemical reaction networks. Parameter sensitivity analysis is a valuable tool that can be used to study robustness properties, for drug targeting, and many other purposes. However its application to stochastic models has been limited when Monte Carlo methods are required due to extremely high computational costs. SPSens provides efficient, state of the art sensitivity analysis algorithms in a single software package so that sensitivity analysis can be easily performed on stochastic models of biochemical reaction networks. SPSens implements the algorithms in C and estimates sensitivities with respect to both infinitesimal and finite perturbations to system parameters, in many cases reducing variance by orders of magnitude compared to basic methods. Included among the features of SPSens are serial and parallel command line versions, an interface with Matlab, and several example problems. AVAILABILITY: SPSens is distributed freely under GPL version 3 and can be downloaded from http://sourceforge.net/projects/spsens/. The software can be run on Linux, Mac OS X and Windows platforms.

A pathwise derivative approach to the computation of parameter sensitivities in discrete stochastic chemical systems.

Characterizing the sensitivity to infinitesimally small perturbations in parameters is a powerful tool for the analysis, modeling, and design of chemical reaction networks. Sensitivity analysis of networks modeled using stochastic chemical kinetics, in which a probabilistic description is used to characterize the inherent randomness of the system, is commonly performed using Monte Carlo methods. Monte Carlo methods require large numbers of stochastic simulations in order to generate accurate statistics, which is usually computationally demanding or in some cases altogether impractical due to the overwhelming computational cost. In this work, we address this problem by presenting the regularized pathwise derivative method for efficient sensitivity analysis. By considering a regularized sensitivity problem and using the random time change description for Markov processes, we are able to construct a sensitivity estimator based on pathwise differentiation (also known as infinitesimal perturbation analysis) that is valid for many problems in stochastic chemical kinetics. The theoretical justification for the method is discussed, and a numerical algorithm is provided to permit straightforward implementation of the method. We show using numerical examples that the new regularized pathwise derivative method (1) is able to accurately estimate the sensitivities for many realistic problems and path functionals, and (2) in many cases outperforms alternative sensitivity methods, including the Girsanov likelihood ratio estimator and common reaction path finite difference method. In fact, we observe that the variance reduction using the regularized pathwise derivative method can be as large as ten orders of magnitude in certain cases, permitting much more efficient sensitivity analysis than is possible using other methods.

Integral tau methods for stiff stochastic chemical systems.

Tau leaping methods enable efficient simulation of discrete stochastic chemical systems. Stiff stochastic systems are particularly challenging since implicit methods, which are good for stiffness, result in noninteger states. The occurrence of negative states is also a common problem in tau leaping. In this paper, we introduce the implicit Minkowski-Weyl tau (IMW-τ) methods. Two updating schemes of the IMW-τ methods are presented: implicit Minkowski-Weyl sequential (IMW-S) and implicit Minkowski-Weyl parallel (IMW-P). The main desirable feature of these methods is that they are designed for stiff stochastic systems with molecular copy numbers ranging from small to large and that they produce integer states without rounding. This is accomplished by the use of a split step where the first part is implicit and computes the mean update while the second part is explicit and generates a random update with the mean computed in the first part. We illustrate the IMW-S and IMW-P methods by some numerical examples, and compare them with existing tau methods. For most cases, the IMW-S and IMW-P methods perform favorably.

Efficient computation of parameter sensitivities of discrete stochastic chemical reaction networks.

Parametric sensitivity of biochemical networks is an indispensable tool for studying system robustness properties, estimating network parameters, and identifying targets for drug therapy. For discrete stochastic representations of biochemical networks where Monte Carlo methods are commonly used, sensitivity analysis can be particularly challenging, as accurate finite difference computations of sensitivity require a large number of simulations for both nominal and perturbed values of the parameters. In this paper we introduce the common random number (CRN) method in conjunction with Gillespie's stochastic simulation algorithm, which exploits positive correlations obtained by using CRNs for nominal and perturbed parameters. We also propose a new method called the common reaction path (CRP) method, which uses CRNs together with the random time change representation of discrete state Markov processes due to Kurtz to estimate the sensitivity via a finite difference approximation applied to coupled reaction paths that emerge naturally in this representation. While both methods reduce the variance of the estimator significantly compared to independent random number finite difference implementations, numerical evidence suggests that the CRP method achieves a greater variance reduction. We also provide some theoretical basis for the superior performance of CRP. The improved accuracy of these methods allows for much more efficient sensitivity estimation. In two example systems reported in this work, speedup factors greater than 300 and 10,000 are demonstrated.

The numerical stability of leaping methods for stochastic simulation of chemically reacting systems.

Tau-leaping methods have recently been proposed for the acceleration of discrete stochastic simulation of chemically reacting systems. This paper considers the numerical stability of these methods. The concept of stochastic absolute stability is defined, discussed, and applied to the following leaping methods: the explicit tau, implicit tau, and trapezoidal tau.